Abstract 



Davydov's model of solitons in a-helix protein chains is shown to display features of 
self-organized criticality (SOC), i.e., power law behaviour of correlations in space and 
1/f noise, as a consequence of considering random peptide group displacements from 
their (periodic) equilibrium positions along a chain. This may shed light on a basic 
mechanism leading to obtain flicker noise in a-helix protein chains and to predict a SOC 
regime in biomolecular structures from first principles. We believe our treatment of 1/f 
noise to be of some relevance to recent findings due to Voss on DNA [Phys. Rev. Lett. 
68, 3805 (1992)]. 
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The concept of solitons found a fascinating novel application into biological phenom- 
ena since Davydov |L| introduced a cubic non-linear Schrodinger equation to describe 
energy transport in molecular chains such as a-helices. The topological and dynamical 
stability of Davydov solitons -i.e., conservation of form and velocity after interaction, 
respectively- is related to the spontaneous (local) symmetry breaking of the protein 
molecules. A dynamical balance between the dispersion due to the resonance interac- 
tion of intrapeptide dipole vibrations and the nonlinearity of the interaction of such 
vibrations with the (time-dependent) local displacements from the equilibrium positions 
of the peptide groups -say, (3- is believed to play a crucial role for the transport of energy 
released in the hydrolisis of the adenosine triphosphate (ATP) molecule 0. 

Because of this, Davydov solitons -as well as alternative soliton models to account for 
nonlinearity in quasi ID bio molecular chains (see, e.g., 0, §)- are an active research 
field. Recent progress has been directed to relate excitations in Davydov systems to 
Bose-Frohlich condensation phenomena || |7j and to include more appropriately quan- 
tum and temperature effects on the lifetime of the soliton |8|, [| . Even though a large 
amount of theoretical work has been accumulated in the past two decades following these 
Davydov's pioneering ideas, there is still a number of open problems. Firstly, there is as 
yet no experimental verification of Davydov solitons [K| PH . Indeed Davydov's model 
might look rather unrealistic || , but from the standpoint of physics it is not unreason- 
able to study it as an effective way towards achieving better understanding of complex 
molecular systems. Secondly, to our knowledge, the case of having random peptide group 
displacements (3 due to local disorder has not been fully considered. 

To this end the recent proposal by Voss |12] (see also |13|, |14[) that long-range cor- 



relations and 1/f noise can be detected in DNA sequences when viewed as random 



processes, gives one inspiration for carrying out an investigation on random (3. The idea 
that nucleotide bases in strands of DNA may be correlated over several thousands of 
base-positions, opens the way to explore how, and to which extent, solitons might be 
useful to gain physical insight into Voss's relevant finding on DNA. It is tempting to 
analyse this phenomenon from a first principle formalism based on soliton physics and, 
from this, to understand that class of non-linear dynamical systems which drive them- 
selves into a statistically stationary critical state, the so-called self-organized criticality 
(SOC) |15| ) with no intrinsic length or time scale, where the systems exhibit power law 



(fractal) behaviour and generate flicker noise [T(| . 

In this work we shall make a small step towards such a theory to show that Davydov 
solitons in quasi-lD a-helix chains at K, if they exist, might display features of SOC 
as a consequence of assuming random peptide group displacements from their (periodic) 
equilibrium positions along a chain. This suggests a possible physical mechanism to 
understand SOC from an ab-initio basis in terms of interactions. In addition our work 
complement recent ideas put forward by one of us [|T7]], regarding an analytical contin- 



uous probability theory for SOC, in that we propose here that random peptide group 
displacements might be one possible mechanism to generate 1/f noise in macromolecular 
chains. 

We start considering the simplest Davydov's protein molecule, i.e. a = 1 0, which 
we briefly describe below for completeness. The basic idea is that the amide-/ vibrational 
energy is coupled through interaction with acoustic phonons. The molecular chain has 
N (3> 1) molecular units of mass M (~ 114 a.m.u.) placed at positions: x n = nR + u n , 
where n is an index that counts unit cells (in the hydrogen-bonded direction), m„(-C R) 
are small displacements from the equilibrium positions nR caused by internal molecular 



motion. The molecular groups are taken as the peptide subunits of a-helix protein 
polymers and R as the length of the amide's hydrogen bond (~ 3A). 

The Hamiltonian operator H for the collective degrees of freedom resulting from the 
interaction of intramolecular amide- J modes, (C = O stretching), and the lattice motion 
of such a chain is 

H = H vib + H ph + Hi n t ■ (1) 
In the above, 

Hyib = ^2{e B n B n — J(B n B n+ i + B} l+1 B n )} , (2) 

n 

where B n and B n are boson creation and annihilation operators for the vibrational 
excitation at the nth site associated with the amide-/ dipolar oscillator -having the 
quantum energy eo ~ 0.205 eV. The dipolar resonant interaction is only considered 
between the nearest neighbor molecules, i.e., J = 2d 2 / B? (= 9.67 x 10~ 4 eV from 
infrared spectra), where d (= 0.29 D) is the dipolar electric moment aligned along the 
C = O bond. 

The second term of Eq. (JJ) describes the longitudinal harmonic oscillations (phonons) 
of the chain, which in second quantized form is written as a sum over the momentum 
normal modes q, namely 

H ph = J2^ q (b% + b > ( 3 ) 
q A 

where 6j and b q are phonon creation and annihilation operators. Q q is given by the 
dispersion equation VL 2 q = 4(z^/i? 2 ) sin 2 (^qR) with v a = R^Jw/M and w an elasticity 
coefficient (~ 76 N/m). 
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The last term in the collective Davydov Hamiltonian is the nonlinear interaction 
between the vibrational dof's and the phonon dof's, namely 

H im = ^tX x(q)e l * nR BlB n (b g + bl q ) , (4) 

V^V q,n 

where x(q) = X*(~Q) = ^x(jMr) 1 ^ 2 sin(gi?). The physical meaning of the nonlinear 
coupling parameter \ is to characterize change effects of the amide-/ bond energy per 
some unit extension of the hydrogen bond. This parameter is expressed as \ = deo/dR, 
whose numerical value lies in the interval: x = 3 — 6.2 x 1CT 11 N. 

For the wavefunction of the above Schrodinger equation, Davydov wrote 

n^(t)»=E a »(*)«p( ff (*)) s iii°» > ( 5 ) 

n 

where ||0)) is a generalized vacuum state in the collective space of vibrational-phonon 
dof's and a n (t) are normalized functions. The quantum-mechanical phase a(t) is written 
as a(t) = — | Ytqifiqnit)^ — /3* n (t)b g ], which is termed the D 1 ansatz. It is a superposition 



of tensor products of single-exciton states and coherent phonon states |18|. There also 
exists a classical displaced oscillator ansatz (or D2) as well as a modified m — D\ state 

E3- 

For the complex functions a n (t) and the real functions I3 qn {t) of Eq.(||) -which char- 
acterize the vibrational state and the displacement from equilibrium of a single molec- 
ular unit at site n, respectively- Davydov obtained a system of two coupled discrete- 
differential equations expressing their 'classical' Hamiltonian evolution. In the long 
wavelength approximation such discrete functions are replaced by two continuous limits 
at(C) an d /9(C)) where ( = nR — v s t, subject to the constraint that the propagation veloc- 
ity v s is constant (i.e., stationary propagation). A further approximation is to take nR 
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as any point x along the chain axis. Besides (3, the function p(£) = —R d/3 gp is also de- 
fined to characterize the (infinitesimal) difference of displacements of nearest neighbour 
molecules within the chain and modulation of the chain due to vibrational solitons. 

Solutions of the energy transport within Davydov's Hamiltonian treatment in the 
continuous and subsonic regime (i.e., s 2 = w s 2 /^<l) are as follows: 

a(C) = ^e^T^-^cosh- 1 ^) , (6) 

p(C) = (7) 

w(l — s z ) 

where Q = 2 wh!\i-s 2 ) • ^ ne Quantities /x, E s , and s are respectively given by /i = 
X7M1 -s 2 ), E s = e -2J+^- Jfi 2 /3, and s = v s /v a . 

In the above a(C) is essentially the Davydov vibrational soliton, which is a quasi- 
localized structure with size of the order R/fi, propagating with velocity v s and trans- 
fering vibrational energy eo- Furthermore, p(() is the hydrogen bond soliton whereas 
/3(C) is a kink soliton of the displacements of the peptide groups from their (periodic) 
equilibrium positions. During their (constant- v s ) movement, the three solitons (usually 
denoted as Si, R\ and Ki, respectively) strongly influence each other in a way which 
is still an open problem. Of special interest to us here is the K\ solution, which might 
be considered as a spatial continuous representation of D'Alambert type (i.e., of the 
form f3(x — v s t)) of a product of coherent states, that is needed to obtain the nonlinear 
Schrodinger dynamics p0| ). 

Let us include next, to the above Davydov model, a special type of random disorder 
that may result from several sources, such as radiation and others, and that fluctuates 
in time. We shall come back to this latter on. In the ^-coordinate frame, the central 



position of the kink K\ can be fixed as the spatial origin and small time excursions of 
its origin will be allowed. Then it is straightforward to estimate time correlations via 
the (dimensionless) noise power spectrum of the temporal evolution of such fluctuations 



where r — > oo and ((x,t) is the random variable. Figure 1 shows the results obtained for 
the spectral density Sp when using x — > x Q ± Ax, such that we set x Q = and choose Ax 
to vary randomly between ±1 for simplicity. These results have been computed using a 
standard fast fourier transform algorithm. On averaging over thousand ensambles, 2 13 
(unit-)time steps within the interval < t < 10 4 are found to suffice in order to achieve 
less than 5% of error in calculations. As revealed from this figure the noise power 
spectrum of the Davydov (3 displacement, when assumed to be a function of random £, 
shows a clear manifestation of 1/f noise. We believe this treatment of flickler noise to 
be somehow of relevance to the findings due to Voss on DNA (12| . In particular, we can 
superpose a peaked structure in the 1/f power spectrum of (3 by periodically resetting 
Ax to naught as shown in the high-frequency region of Fig.l. 

From these results we can see that the transition to a steady state, is strongly re- 
lated to the peptide group fluctuations from their initial (i.e., t Q = 0, x Q = 0) conditions. 
According to such complex dynamics, the physics behind this phenomenon -i.e., exis- 
tence of time correlations- is due to the random disorder we have introduced and whose 
possible cause will be analysed below. 

Let us first also investigate the space domain by fixing the time scale (i.e., by assum- 




1// by 




(9) 
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ing now static solutions). We consider once more the random variable £ to characterize 
random disorder of the peptide group displacements from their equilibrium positions 
with a density probability <fi proportional to /3(C) W% - Since the uniform probability 
distribution function of having random events can be written as 

<?(Ca) " <?(Ci) = nG < C < C2} w 0(C) dC , (10) 

where {} indicates the function interval. Then, this integral over the limits (2 = ^1(0 > 
C + A2C0 = Ci gives 

0(AiC«O - 0(C + A 2 Co) = / Al? ° /3(C) dC' = -r(C) , (11) 

■/C+A2C0 

with /3 given in Eq.([S|), (s 2 <C 1), and Co being a constant. The minus implies that the 
functions Q are here assumed to satisfy the condition + A 2 Co) > £/(AiC) for C 7^ 0, 
which do not need to be defined, whereas the free parameters A» (z=l,2) will restrict the 
range of £. 

The afore- mentioned integration limits lead to the condition A 2 — Ai + C < as 



discussed in our previous work |17| . If A 2 7^ Ai, then we get (in terms of r(0)) 



r(C*) « (1 + -^— ){r(0) + C o lnco^A 2 } 
A 2 — Ai 

-Co{ - Incog/iA! + lncog/i(A 2 + C*)} • (12) 
A 2 — Ai 

where C* = C/C, r(0) = Co(A 2 - A x )(l + T A ) and T A = (IncoshXt - In cos/iA 2 )/(A 2 - 

In Fig. 2 we show the dependence of the normalized probability distribution function 
r on the reduced variable (* for values of Ai = 3, A 2 = —8 and r(0) = 1 which, in 
turn, determine the value of ( a . In our calculations, i.e. using Eqs.(g) and Eq. flT2|) , we 
have also reduced W ^L S 2^ ^ 1 an d Q( oc 1 for simplicity. The choice of r(0) allows 
to normalize r(C*) and to mimic features of SOC, namely a power-law behaviour in 



space correlations (within the range < < 8), provided (* is associated with the 
Zog-function of a measured random event. In fact, we have that t((* — > 0) — ► 1 and 
r(Ai - A 2 ) = 0. Whereas if A 2 + C* « then In cosh(X 2 + (*) « 0, hence r of Eq. ([12]) 
depends linearly on (* for values C* < — A 2 . 

To see more clearly possible power-law features in the behaviour of r over an extensive 
range of values of (* < 8, we calculate the linear derivative of r(£*) which is also plotted 
in Fig. 2 by dotted lines. For the smallest displayed values of (* this function converges 
to a constant negative value, revealing in this way the constant nature of the negative 
slopes in the (t-£*) curves. In view of these features of the derivative of r, the second 
derivative has also been included in Fig. 2. It presents a sharp peak around the inflection 
point of J^r, thus indicating the range of validity of such a power-law behaviour. The 
cutoff in the £*-axis for the r curve is related to the system size or integration limit 
(i.e., Ai — A 2 ~ 11). Accordingly, it reflects the range of long-range correlations in 
the space domain. In view of these results we can interpret f3 as a nonequilibrium 
order parameter of a transition from the D\ phase of the chain to the Frohlich phase 
| 2T| . That is, a phase transition from the dynamical balance between the intramolecular 



excitations (and their exchange) and the longitudinal excitations of the linear Davydov 
chain, to the Frohlich nonthermal excitations of longitudinal polarization modes arising 
in a far-from-equilibrium regime supported by the flow of energy. 

We focus now on some possible effects that have been usually neglected when taking a 
continuous limit within the framework of the Davydov theory but which, in our opinion, 
may become important as being the source for generating randomness. These effects 
include radiation, discreteness-chaos, and disorder. Radiation effects have been first 
discovered in the context of lattice topological solitons (dislocations), and subsequently 
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proved to exist also in the case of dynamical solitons |22]| . In fact, subsonic kinks in 
a monoatomic chain permanently radiate small amplitude oscillations. Besides this, 
subsonic kinks (as well as supersonic ones) in diatomic chains lose energy in this process 



] . Chaotic effects, due to lattice discreteness have been discussed in Ref . [f23|] within a 
vibron model, with on-site potentials, which is very similar to the Davydov model. Such 
an effect can be understood as a perturbation of an integrable system. It is a chaotic 
effect that might imply a random spatial arrangement of stationary solitons. The effect 



of disorder, on the other hand, is a very important feature of many nonlinear chains |24 
Considering the 20 different natural amino acids in real peptide chains, each of different 
mass, one may think of a sequence of random masses to have important effects on 
soliton propagation. This has been observed by Forner |25| by investigating sequences 
of masses, spring constants, nonlinear coupling constants, heat bath, and disorder in 
the dipole coupling. Another crucial effect to consider when applying the continuous 
limit is the impurity disorder in which various types of kink-impurity interactions may 



be possible p6 |. We believe all of these phenonomena to generate (a sort of intrinsic) 



randomness for the displacements (3 along a single molecular chain {i.e., a = 1), which 
leads to obtain SOC features as we have discussed in this work. We add that the present 
ideas may be easily extended to chains with a > 1, since the soliton solutions of the 
Davydov Hamiltonian -i.e., S±, R±, K\- keep their form unchanged [|l|]. 

We have thus tried to interpret, within the simplest Davydov soliton theory, ran- 
domness by combining features of both discreteness and disorder effects. These consid- 



erations may be useful to understand long-range correlations in biological systems f27 
and flicker noise, possibly including DNA [|l!J. We have predicted a SOC dynamical 
regime in biomolecular systems from first principles. Indeed we have been able to derive 



such a SOC regime as a consequence of random (3 which has been interpreted in terms 
of a nonequilibrium phase transition between the D\ soliton state of the chain and a 
probable Frohlich- condensed phase. The present SOC regime might be seen as a proof 
of the soliton stability against random disorder. In turn, this may be an indication of 
the nonequilibrium nucleation of Frohlich domains along a single chain. Of course, to 
characterize randomness of the peptide group displacements with a density probability 
oc j3(() may be seen heuristic. But, in view of the results obtained, this can be con- 
sidered as one of the simplest reasonable ways to relate local properties to macroscopic 
behaviour. To this end, note that a relation of this kind has also been used within the 
context of pattern formation |28fl . 

The authors thank Dr. J. Langowski (EMBL-Grenoble) for a valuable discussion and 
Prof. Yu Lu (ICTP-Trieste) for encouragement. 
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Figure captions 

• Fig.l: Noise power spectrum for the temporal evolution of the random fluctuations 
in (3(C) including periodical resetting. 

• Fig. 2: (Full line) probability distribution function r of having random events vs. 
reduced variable (* using Ai > A 2 , such that A 2 < 0, and r(0) = 1. (Dotted lines) 
first and (smallest dotted lines) second derivative of r with respect to (*. 
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